Complexity in dislocation dynamics: model 
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We propose a numerical model to study the viscoplastic deformation of ice single crystals. 
We consider long-range elastic interactions among dislocations, the possibility of mutual 
annihilation, and a multiplication mechanism representing the activation of Frank-Read 
sources due to dislocation pinning. The overdamped equations of motion for a collection 
of dislocations are integrated numerically using different externally applied stresses. Using 
this approach we analyze the avalanche-like rearrangements of dislocations during the dy- 
namic evolution. We observe a power law distribution of avalanche sizes which we compare 
with acoustic emission experiments in ice single crystals under creep deformation. We em- 
phasize the connections of our model with non-equilibrium phase transitions and critical 
phenomena. 



I. INTRODUCTION 

The viscoplastic deformation of crystalline mate- 
rials, such as ice single crystals, involves the motion 
of a large number of dislocations. Dislocations may 
be incorporated into a crystal in its growth process, 
or under deformation conditions. They can pene- 
trate into the material from the sample surfaces, or 
be generated by various mechanisms, as for exam- 
ple, in what is usually called a Frank-Read source, 
activated by the pinning of a dislocation loop 1 . 
The elastic force between a pair of dislocation lines 
decreases algebraically with the interline distance; 
it can be attractive or repulsive, depending on the 
orientation of their respective Burgers vectors; and 
eventually a pair of dislocations may annihilate at 
relatively short distances. As a consequence of all 
of these features, the collective behavior of a large 
number of these defects appears to be an amazingly 
rich problem, suitable to be studied from very dif- 
ferent points of view. 

A new insight into the classical issue of vis- 
coplasticity was recently opened by Weiss and 
Grasso after their experimental study of the 

acoustic activity of ice crystals under creep. The 
complex character of the collective dislocation dy- 
namics reveals itself in experiments of acoustic 
emission (AE), as reported in our companion pa- 
per 1^ in this volume. As a matter of fact, when a 
material is deformed under constant load (creep ex- 
periment) and dislocation motion is the dominant 
mechanism for viscoplastic deformation, a constant 
strain-rate regime usually follows after the initial 
transient stage. Orowan's relation 7 — pmbv is 
known to prevail under such conditions, where 7 
is the strain of the sample, pm is the density of mo- 
bile dislocations, b is the Burgers' vector, and v is 
the mean velocity of the dislocations. Obviously, 
this is a mean-field relation which neglects tempo- 



ral and spatial fluctuations of both the density and 
the velocity fields. As a result of their interactions, 
however, dislocations tend to move cooperatively 
giving rise to a rather complex and heterogeneous 
slip process. Dislocations move in groups to form 
slip bands. Moving dislocations can pile-up against 
stable dislocation configurations such as walls or 
boundaries, which may eventually break apart. 

Given the amplitude threshold and the frequency 
range accessible to the experimental apparatus, the 
AE signals detected seem to correspond to the syn- 
chronous motion of several dislocations, likely to oc- 
cur during the breakaway of a pile of these defects, 
or the activation of a multiplication source. Vari- 
ous measurements of the acoustic activity recorded 
during a stress-constant step show that the AE sig- 
nal takes place in the form of bursts distributed ac- 
cording to a power law. This behavior is certainly a 
consequence of the collective motion of dislocations 
which spontaneously gives rise to an avalanche-like 
dynamics, typical of slowly driven dissipative sys- 
tems. The power law distribution provides evi- 
dence of scale-free cooperative behavior whose ori- 
gin could be ascribed to nonequilibrium continuous 
phase transitions ||^ or self-organized criticality . 

The AE experiments, however, have only access 
to information resulting from the interplay of vari- 
ous magnitudes. Thus, the physical interpretation 
of the generated AE waves remains a major diffi- 
culty, and constitutes the main motivation of this 
work. 



II. DESCRIPTION OF THE MODEL 

Our goal is to characterize the viscoplastic de- 
formation of ice or a similar crystalline material 
from the perspective of nonequilibrium statistical 
mechanics. To do so, we propose a simplified model 
to study the dynamics of a collection of dislocations. 



Ice crystals deform essentially by slip on the basal 
plane (0001) [xy plane), i.e. the motion of disloca- 
tion lines or loops takes place by gliding on the xy 
planes. The simplest Burgers vectors b of a disloca- 
tion in a crystal of hexagonal ice are the primitive 
cell vectors of the conventional hexagonal lattice. 
We start studying a two-dimensional (2d) model 
representing a cross section of the crystal which is 
perpendicular to the basal planes and parallel to 
one of these lattice vectors, that is, for example, the 
xz plane. In this way, the dislocations constrained 
to move in this plane have all Burgers vectors par- 
allel to the chosen lattice vector b — {b, 0, 0) and 
move along fixed lines parallel to the x axis. We 
also consider that all N dislocations are of edge 
type, and that, on average over many realizations, 
the number of dislocations with positive and nega- 
tive Burgers vectors is the same. 

Several 2d-models containing similar basic ingre- 
dients have been proposed in the literature in the 
last few years ]^-p^. A basic feature common to 
most models, is that dislocations interact with each 
other through the long-range elastic stress field they 
produce in the host material. An edge dislocation 
with Burgers vector b located at the origin gives rise 
to a shear stress cr* at a point {x, z) of the form 



where D = /i/27r(l — cr) is a coefficient involving 
the shear modulus /i and the Poisson's ration a for 
the material. In our model, we further assume that 
Peierls stress is zero and that the dislocation ve- 
locities are linearly proportional to the local stress. 
Experimental evidence supports such a relationship 
for low stress conditions, which is indeed the case 
in our model. Accordingly, the velocity of the nth 
dislocation, if an external shear stress tr'^ is also ap- 
plied, is given by 

v^^b^{Y,<^~^% (2) 

As the number of dislocations in any real crystal 
exceeds by far the number of defects we can han- 
dle in a computer, one usually introduces periodic 
boundary conditions (PBC) to effectively extend 
the size of our system. To avoid the discontinuities 
arising from truncating long-range elastic interac- 
tions (see Eq. |l|) , we resort to the Ewald summation 
method. In this way, we have exactly accounted 
for the interaction of a dislocation with all the in- 
finite periodic replicas of another dislocation in a 
finite cell of dimensions L x L. Contrary to what 
is stated in Ref. we do not find any spurious 
results coming from the implementation of PBC in 
this fashion. 

When the distance between two dislocations is of 
the order of a few Burgers vectors, the high stress 



and strain conditions close to the dislocation core 
invalidate the results obtained from a linear elastic- 
ity theory (i.e. Eq .(|l|)). In these instances, phe- 
nomenological nonlinear reactions describe more ac- 
curately the real behavior of dislocations in a crys- 
tal. In particular in our model, we account for the 
annihilation of dislocations with opposite Burgers 
vectors when the distance between them is shorter 
than 26. Thus the core of one dislocation in our 
model has a radius of size b. 

Another important feature of any computer 
model is the implementation of a mechanism for the 
multiplication of dislocations in the sample. It is 
widely believed that the Frank-Read mechanism [Q 
is the most relevant for a gliding process of dis- 
locations under creep deformation. Indeed Frank- 
Read sources (FRS) have been observed in ice. In a 
FRS multiplication occurs by pinning of a disloca- 
tion segment on the basal planes due, for example, 
to a defect in the crystal, or to the very presence 
of dislocation dipoles, piles, and walls. If the lo- 
cal stress concentration is less than a critical value, 
the pinned segment bows out by glide. Beyond this 
critical value, the dislocation segment wraps around 
itself, creating a new dislocation loop and restoring 
the original configuration. Thus a sequence of loops 
forms continuously from the source until the local 
shear stress drops below the activation value. In or- 
der to simulate this mechanism, we split our system 
of size L into smaller cells where, at each time step, 
we keep track of both the local fraction of pinned 
dislocations and the stress. Pairs of dislocations of 
opposite Burgers vectors are generated in a given 
cell if (i) there are pinned particles and (ii) if the 
local stress is large compared to a threshold value 
(the pair creation probability is proportional to the 
local stress). Each pair of dislocations reduces the 
local value of the stress, and we keep on generating 
pairs until the local stress drops below the threshold 
value. In such a way, we simulate the procession of 
dislocations often seen experimentally when a FRS 
is activated in a material. 

To follow the time evolution of N dislocations, 
we integrate numerically the N coupled equations 
dxn/dt^b„{Y,^_^^ <m - <^''), where ct^„ is given 
by Eq. m, using an adaptive step size fifth- order 
Runge-Kutta algorithm. Simulations start from 
a configuration of N point-dislocations randomly 
placed on a square cell of size L. So far, we have 
considered two different box sizes L = 2006, and 
L = 3006, with an initial number of dislocations 
iVo = 800, and iVo = 1500, respectively. We let the 
system relax in the absence of external stress, until 
it reaches a metastable arrangement. Annihilation 
and multiplication processes imply that the num- 
ber of dislocations N in our system is not fixed in 
the course of time. After the system has reached a 
metastable arrangement, the volume fraction of dis- 
locations 4) = Nnb'^ / (Lb)'^ ranges between 1 — 5%. 



2 



Once in these conditions, we apply a small external 
shear stress and keep track of the various quantities 
describing the dynamics of the dislocations. Aver- 
ages are made over at least 50 realizations starting 
with different random initial conditions. 



III. CREEP DYNAMICS 

In the evolution of our model, dislocations build 
up complex and highly fluctuating patterns contain- 
ing, for example, dislocation dipoles and walls, sim- 
ilar to those experimentally observed ]T3| , |l^ . 

If the external stress applied is not high enough 
to activate a FRS our system simply relaxes 
very slowly. After a short initial transient, the 
root-mean-square velocity vfY^'^ and the mean 
strain rate d'y/dt ~ ^i^iVi enter into a regime of 
power-law relaxation. In Fig. |^ we show in a dou- 
ble logarithmic scale the relaxation of d'y/dt. The 
power-law exponent is around —0.65. This scale- 
free behavior seems to correspond to what in the 
literature is known as Andrade's power law creep, 
where the strain rate decays as t"^/"^, and, conse- 
quently, the strain grows like t^/"^. The inset of 
Fig. |l| shows the strain curve. Our best fit is ob- 
tained with an exponent of 0.28. 




FIG. 1. Mean strain rate as a function of time in 
our model of collective dislocation dynamics. The curve 
shown in the inset represents the global strain. 

One can further investigate the behavior of other 
quantities involved in the collective dislocation dy- 
namics, such as the moving dislocation density, the 
annihilation rate, as well as the trajectories of the 
dislocations. These results will be presented in de- 
tail elsewhere [^2). 

If we increase the external applied stress, the sys- 
tem starts to generate dislocation pairs through the 
mechanism described in the previous section. At 
this point, the system crosses over to a regime in 
which the strain rate relaxes to a constant value. 



i.e. the linear creep regime. In this regime, we can 
monitor observables such as the average number of 
dislocations and the mean velocity V = V' \vi\ of 
all moving dislocations Nm- In Figure 0, we see 
that V exhibits an intermittent and bursty behav- 
ior that signals the presence of collective dislocation 
rearrangements. Noticeably, the highest activity is 
not necessarily in a one-to-one correspondence with 
the total number of dislocations N present in the 
system. 
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FIG. 2. Number of dislocations A'^ and mean velocity 
y as a function of time in a given run of our model with 
generation of dislocation pairs in FRS's. 

The features of this signal resemble those ob- 
served in apparently very different systems , such as 
fracturing of wood and concrete |I6|| , Barkhauscn 
effect [|l7|l , and flux lines in high-Tc superconduc- 
tors [ |I8| , and induces us to think in terms of 
avalanches. To define what an avalanche is in our 
system, we introduce an offset value of the signal to 
eliminate what we call the background noise. The 
latter can be thought of as the lower experimen- 
tal resolution, but at the same time corresponds to 
the incoherent motion of Nm independent disloca- 
tions. Thus, when the signal exceeds this offset, 
we will have an avalanche. The avalanche dura- 
tion T and magnitude s = ^ V increase until the 
signal drops down below the offset value. The cu- 
mulative distribution Pc{s) of the avalanches ob- 
tained after averaging over several realizations is 
depicted in Figure |^ for the two system sizes stud- 
ied L — 200, 300. We recover a clear power law 
distribution {Pc{s) ^ s^^) extending over close to 
two decades. The exponent r ~ 0.6 is in reasonable 
agreement with experimental data The univer- 
sal properties of the distribution does not depend on 
the particular choice of the offset value. The distri- 
bution cut-off for large values of s is due to the finite 
size of the sample. (As one would expect, the cut- 
off is scaling accordingly to the system size.) This 
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clearly points out the presence of a very large (or 
infinite) characteristic size for the acoustic events. 
It is worth remarking that larger stresses intro- 
duce a characteristic scale in the process. A more 
detailed study as a function of the applied stress is 
in progress 




FIG. 3. The cumulative avalanche size distribution 
Pc for two system sizes L — 200, 300 represented in a 
double logarithmic scale. The curves in the inset repre- 
sent the raw distributions P{s). 

IV. CONCLUSIONS AND DISCUSSION 

The numerical investigation of the present model 
provides striking evidence for the collective criti- 
cal behavior of dislocation motion under external 
stress. We have found two different regimes de- 
pending on the externally applied stress: i) For low 
stresses, the mean dislocation velocity and the mean 
strain-rate relax according to Andrade's power law 
creep, ii) For higher stress values, the system acti- 
vates FRS's to which it responds exhibiting singu- 
lar behavior in the guise of avalanches distributed 
over many length scales. The avalanche distribu- 
tion is a power law which signals the absence of 
any characteristic size in the process. Avalanche 
dynamics is the rule rather the exception in slowly 
driven disordered systems . Under the exter- 

nal drive (the stress in the present case), the sys- 
tem jumps between metastable or pinned configura- 
tions in which the dynamics is virtually frozen. In 
the limit of a very slow driving the disordered en- 
ergy landscape is explored quasistatically and the 
response function exhibits critical properties pof . 
Typically, a basic ingredient for this behavior is the 
presence of quenched disorder acting as the source 
of pinning in the system. Noticeably, the system 
under study does not contain any external source of 
disorder. Pinned states are due to the various struc- 
tures such as dipoles, piles, and dislocation walls, 
that play the role of self-generated pinning centers 



that create the pinning force landscape. The new 
scenario poses many new and interesting questions 
for a definitive identification and understanding of 
the critical nature of dislocation dynamics. 
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